/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | www.openfoam.com
     \\/     M anipulation  |
-------------------------------------------------------------------------------
    Copyright (C) 2024 AUTHOR,AFFILIATION
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.

Application
    checkPointsExist

Description

\*---------------------------------------------------------------------------*/

#include "fvCFD.H"
#include <chrono> // 函数运行计时
#include "meshSearch.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

// 从csv中读取数据，获得所有网格点
std::vector<std::vector<std::string>> readCSV(const std::string& filename) {
    std::ifstream file(filename);
    std::vector<std::vector<std::string>> data;
    std::string line;

    if (file.is_open()) {
        while (std::getline(file, line)) {
            std::stringstream ss(line);
            std::string item;
            std::vector<std::string> row;

            while (std::getline(ss, item, ',')) {
                row.push_back(item);
            }
            data.push_back(row); // loop through each row
        }
        file.close();
    } else {
        std::cerr << "Unable to open file" << std::endl;
    }
    return data;
}

// 检查输入点是否能够从网格中找到cell, method 1
bool checkPointsInCell(const double (&point)[3], const fvMesh& mesh)
// bool checkPointsInCell(const double (&point)[3])
{
    vector coord(point[0], point[1], point[2]);
    Info << coord << endl;
    // 使用 findCell 函数查找点所在的单元
    label cellId = mesh.findCell(coord);
    return cellId != -1;
}

// 检查输入点是否能够从网格中找到cell, method 2
bool checkPointsInCell_m2(const double (&point)[3], const meshSearch& ms)
// bool checkPointsInCell(const double (&point)[3])
{
    vector coord(point[0], point[1], point[2]);
    Info << coord << endl;
    // 使用 findCell 函数查找点所在的单元
    label cellId_1 = ms.findCell(coord, -1, true);
    label cellId_2 = ms.findNearestCell(coord, 0, true);
    return cellId_1 != -1 || cellId_2 != -1;
}

int main(int argc, char *argv[])
{
    #include "setRootCase.H"
    #include "createTime.H"
    #include "createMesh.H"

    std::vector<std::vector<std::string>> data = readCSV("data.csv");

    // 提取第2、3、4列的值
    std::vector<std::vector<std::string>> result;
    //从第二行开始，忽略第一行
    std::vector<std::string> col;

    for (size_t i = 1; i < data.size(); i++) { 
        col.push_back(data[i][1]); // 第2列
        col.push_back(data[i][2]); // 第3列
        col.push_back(data[i][3]); // 第4列
        result.push_back(col);
    }

    // for(const auto& item : col){
    //     std::cout << item << " ";
    // }

    const size_t num_points = col.size() / 3;
    std::cout << "Num of points: " << num_points << std::endl; 

    double location[3] = {0.0};

    auto start = std::chrono::high_resolution_clock::now();

    //- find cell using method 1
    for(size_t i = 0; i < col.size(); i = i + 3)
    {
        location[0] = std::stod(col[0 + i]);
        location[1] = std::stod(col[1 + i]);
        location[2] = std::stod(col[2 + i]);

        // std::cout << "location[0]: " << location[0] << std::endl;
        // std::cout << "location[1]: " << location[1] << std::endl;
        // std::cout << "location[2]: " << location[2] << std::endl;
        std::cout << "index: " << i << std::endl;

        // 开始计时
        auto start_inloop = std::chrono::high_resolution_clock::now();
        if (checkPointsInCell(location,mesh) == true)
        {
            Info << "Point: " << location[0] << ", " 
                              << location[1] << ", " 
                              << location[2] << "find cell in mesh" << endl;
        }else
        {
            Info << "WARN------Point: " << location[0] << ", " 
                                        << location[1] << ", " 
                                        << location[2] << " not find cell in mesh" << endl;
        }
        // 结束计时并计算运行时间
        auto end_inloop = std::chrono::high_resolution_clock::now();
        // auto duration_inloop = std::chrono::duration_cast<std::chrono::microseconds>(end_inloop - start_inloop).count();
        auto duration_inloop = std::chrono::duration_cast<std::chrono::duration<double>>(end_inloop - start_inloop).count();
        std::cout << "\n Function took " << duration_inloop << " seconds" << std::endl;
    }

    auto end = std::chrono::high_resolution_clock::now();
    auto duration = std::chrono::duration_cast<std::chrono::duration<double>>(end - start).count();
    std::cout << "\n Method1 took " << duration << " seconds" << std::endl;

    // find cell using method 2
    const meshSearch ms(mesh);
    start = std::chrono::high_resolution_clock::now();
    for(size_t i = 0; i < col.size(); i = i + 3)
    {
        location[0] = std::stod(col[0 + i]);
        location[1] = std::stod(col[1 + i]);
        location[2] = std::stod(col[2 + i]);

        // std::cout << "location[0]: " << location[0] << std::endl;
        // std::cout << "location[1]: " << location[1] << std::endl;
        // std::cout << "location[2]: " << location[2] << std::endl;
        std::cout << "index: " << i << std::endl;

        // 开始计时
        auto start_inloop = std::chrono::high_resolution_clock::now();
        if (checkPointsInCell_m2(location,ms) == true)
        {
            Info << "Point: " << location[0] << ", " 
                              << location[1] << ", " 
                              << location[2] << "find cell in mesh" << endl;
        }else
        {
            Info << "WARN------Point: " << location[0] << ", " 
                                        << location[1] << ", " 
                                        << location[2] << " not find cell in mesh" << endl;
        }
        // 结束计时并计算运行时间
        auto end_inloop = std::chrono::high_resolution_clock::now();
        auto duration_inloop = std::chrono::duration_cast<std::chrono::duration<double>>(end_inloop - start_inloop).count();
        std::cout << "\n Function took " << duration_inloop << " seconds" << std::endl;
    }

    end = std::chrono::high_resolution_clock::now();
    duration = std::chrono::duration_cast<std::chrono::duration<double>>(end - start).count();
    std::cout << "\n Method2 took " << duration << " seconds" << std::endl;
    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //



    Info<< nl;
    // runTime.printExecutionTime(Info);

    Info<< "End\n" << endl;

    return 0;
}


// ************************************************************************* //
